!------------------------------------------------------------------------------
! TJU/Department of Mechanics, Fluid Mechanics, Code START
!------------------------------------------------------------------------------
!> @Author: Jianxin Liu shookware@tju.edu.cn
!> @Date: 2024-08-17 20:51:30
!> @LastEditTime: 2024-08-17 20:51:32
!> @LastEditors: Jianxin Liu shookware@tju.edu.cn
!> @Description:
!> @FilePath: /START_NS_3D/src/reconstructor/module_reconstructer_muscl.f90
!> @Copyright (c) 2024 by Jianxin Liu email: shookware@tju.edu.cn, All Rights Reserved.
!------------------------------------------------------------------------------
submodule(mod_reconstructor) submod_reconstructor_muscl

contains

   module subroutine CreateMUSCL(this, length, stencil_width, flg_bounds, limiter)
      class(type_reconstructor_muscl), intent(inout) :: this
      integer, intent(in) :: length
      integer, intent(in) :: stencil_width
      logical, intent(in) :: flg_bounds(2)
      integer, intent(in), optional :: limiter

      integer :: i

      this%length = length
      this%stencil_width = stencil_width
      this%flg_bounds = flg_bounds
      this%scheme_limiter = optval(limiter, LIMITER_NONE)
      if (flg_bounds(1)) then
         this%ista = stencil_width
      else
         this%ista = 0
      end if
      if (flg_bounds(2)) then
         this%iend = length - 1 - stencil_width
      else
         this%iend = length - 1
      end if

      call this%is_p%create(length, [(i, i=0, length - 1)])
      call this%is_m%create(length, [(i, i=1, length - 1), length - 1])

      call this%select_limiter(this%scheme_limiter)

   end subroutine CreateMUSCL

   module subroutine ReconstructMUSCL(this, f, h, s)
      class(type_reconstructor_muscl), intent(in) :: this
      class(type_vector_base), intent(in) :: f
      class(type_vector_base), intent(inout) :: h
      character, intent(in), optional :: s

      class(type_self_vector), pointer :: f_ptr, h_ptr

      call vector_2_selfvector(f, f_ptr)
      call vector_2_selfvector(h, h_ptr)
      call this%reconstruct2(f_ptr, h_ptr, s)

   end subroutine ReconstructMUSCL

   module subroutine ReconstructMUSCL2(this, f, h, s)
      class(type_reconstructor_muscl), intent(in) :: this
      class(type_self_vector), intent(in) :: f
      class(type_self_vector), intent(inout) :: h
      character, intent(in), optional :: s
      real(dp), pointer :: f_array(:), h_array(:)
      real(dp), allocatable :: df_array(:)
      integer :: i, l, ista, iend
      real(dp) :: r, limiter

      !!因为get_self_array_read导致数组从1开始计数
      ista = this%ista + 1
      iend = this%iend + 1
      call f%get_array_read(f_array)
      call h%get_array(h_array)
      allocate (df_array(ista + 1:iend))
      df_array = f_array(ista + 1:iend) - f_array(ista:iend - 1)
      select case (optval(s, "+"))
      case ("+")
         h_array(ista) = f_array(ista)
         do concurrent(i=ista + 1:iend - 1)
            r = df_array(i + 1)/df_array(i)
            limiter = 0.25_dp*((1.0_dp - this%kappa)*this%limiter(r) + (1.0_dp + this%kappa)*this%limiter(1.0_dp/r)*r)
            h_array(i) = f_array(i) + 0.5d0*this%limiter(r)*(f_array(i) - f_array(i - 1))
         end do
         h_array(iend:this%length) = 0.0_dp
      case ("-")
         do concurrent(i=ista:iend - 2)
            r = df_array(i + 2)/df_array(i + 1)
            limiter = 0.25_dp*((1.0_dp - this%kappa)*this%limiter(r)/r + (1.0_dp + this%kappa)*this%limiter(1.0_dp/r))
            h_array(i) = f_array(i + 1) - 0.5d0*this%limiter(r)*(f_array(i + 1) - f_array(i))
         end do
         h_array(iend - 1) = f_array(iend)
         h_array(iend:this%length) = 0.0_dp
      end select
      call f%restore_array_read(f_array)
      call h%restore_array(h_array)
      deallocate (df_array)

      !@NOTE: 应该还有其他的对应的不同limiter的写法，这里是直接参考闫超书上的写法，实际上还有其他的写法，比如CFL3D手册上的函数等等

   end subroutine ReconstructMUSCL2

   module subroutine DestroyUpMUSCL(this)
      class(type_reconstructor_muscl), intent(inout) :: this

      call this%is_p%destroy()
      call this%is_m%destroy()
      this%limiter => null()

   end subroutine DestroyUpMUSCL

end submodule submod_reconstructor_muscl
